%code for calibration of parameters and calculation of mcf by using the model
%without informality
clear all
clc
syms theta
global tauk taun tauc alpha fgdp phi b delta Time

filename = 'data.xlsx';
data = xlsread(filename);

Time=16*330;

[r,c]=size(data);

 for i=1:r
    alpha=data(i,1);
    b=data(i,2);
    delta=data(i,3);
    
    tauk=data(i,4);
    taun=data(i,5);
    tauc=data(i,6);
    
    fgdp=data(i,7);
   
    phi=data(i,9);
    
    
labor=Time/(1+phi+((alpha*(1-tauk)*phi)/((1-alpha)*(1-taun)))*((1-b)/(1-b+delta*b)));   
K1=(((1-tauk)*theta*alpha)/(1/b-1+delta))^(1/(1-alpha));
K=labor*K1;
Y=theta*(K^(alpha))*(labor^(1-alpha));
C=((1-tauk)*alpha*Y+(1-taun)*(1-alpha)*Y-delta*K)/(1+tauc);
leisure=Time-labor;
R=tauc*C+tauk*alpha*Y+taun*(1-alpha)*Y;

%calibration of parameter theta
z0=0.5;
res=fsolve(@(theta) funcs (theta),z0);
%calculation of the values by using calibrated varaibles
leisureres=subs(leisure,{theta},{res});
laborres=subs(labor,{theta},{res});
Kres=subs(K,{theta},{res});
Cres=subs(C,{theta},{res});
Yres=subs(Y,{theta},{res});
Rres=subs(R,{theta},{res});

data(i,11)=eval(leisureres);
data(i,12)=eval(laborres);
data(i,13)=eval(Kres);
data(i,14)=eval(Cres);
data(i,15)=res;
data(i,16)=eval(Yres);
data(i,17)=eval(Rres);


 end
 
 filename='input.xlsx';
 xlswrite(filename,data);
 
 
 clear all
clc


filename = 'input.xlsx';
mcfdata = xlsread(filename);

Time=16*330;

syms tauk taun tauc theta alpha phi delta b
labor=Time/(1+phi+((alpha*(1-tauk)*phi)/((1-alpha)*(1-taun)))*(((1/b)-1+2*delta)/((1/b)-1+delta)));   
K1=(((1-tauk)*theta*alpha)/(1/b-1+delta))^(1/(1-alpha));
K=labor*K1;
Y=theta*(K^(alpha))*(labor^(1-alpha));
C=((1-tauk)*alpha*Y+(1-taun)*(1-alpha)*Y-delta*K)/(1+tauc);
leisure=Time-labor;

[r,c]=size(mcfdata);

R=tauc*C+tauk*alpha*Y+taun*(1-alpha)*Y;
utility=log(C)+phi*log(leisure);
turevRtauc=diff(R,tauc);
turevUtauc=diff(utility,tauc);
turevRtauk=diff(R,tauk);
turevUtauk=diff(utility,tauk);
mcpftauk=-(turevUtauk/turevRtauk);
turevRtaun=diff(R,taun);
turevUtaun=diff(utility,taun);
mcpftaun=-(turevUtaun/turevRtaun);
mcpftauc=-(turevUtauc/turevRtauc);


for i=1:r
  %calculation of marginal cost of public funds
    mcfk=subs(mcpftauk,{alpha,b,delta,tauk,taun,tauc,phi,theta},{mcfdata(i,1),mcfdata(i,2),mcfdata(i,3),mcfdata(i,4),mcfdata(i,5),mcfdata(i,6),mcfdata(i,9),mcfdata(i,15)});
    mcfc=subs(mcpftauc,{alpha,b,delta,tauk,taun,tauc,phi,theta},{mcfdata(i,1),mcfdata(i,2),mcfdata(i,3),mcfdata(i,4),mcfdata(i,5),mcfdata(i,6),mcfdata(i,9),mcfdata(i,15)});
    mcfn=subs(mcpftaun,{alpha,b,delta,tauk,taun,tauc,phi,theta},{mcfdata(i,1),mcfdata(i,2),mcfdata(i,3),mcfdata(i,4),mcfdata(i,5),mcfdata(i,6),mcfdata(i,9),mcfdata(i,15)});
   
    mcfdata(i,18)=eval(mcfk);
    mcfdata(i,19)=eval(mcfc);
    mcfdata(i,20)=eval(mcfn);
 
end

filename = 'output.xlsx';
xlswrite(filename,mcfdata);

 